Spatial Data as Matrices and Rasters

HES 505 Fall 2022: Session 10

Matt Williamson

Today’s Plan

Image Source: USGS

Objectives

  • By the end of today, you should be able to:

    • Describe the raster data model and its representation in R

    • Access the elements that define a raster

    • Build rasters from scratch using matrix operations and terra

Defining the raster data model

Image Source: QGIS User’s manual

The Raster Data Model

Defining the Raster Data Model

  • Vector data describe the “exact” locations of features on a landscape (including a Cartesian landscape)

  • Raster data represent spatially continuous phenomena (NA is possible)

  • Depict the alignment of data on a regular lattice (often a square)

    • Operations mimic those for matrix objects in R

  • Geometry is implicit; the spatial extent and number of rows and columns define the cell size

Types of Raster Data

  • Regular: constant cell size; axes aligned with Easting and Northing

  • Rotated: constant cell size; axes not aligned with Easting and Northing

  • Sheared: constant cell size; axes not parallel

  • Rectilinear: cell size varies along a dimension

  • Curvilinear: cell size and orientation dependent on the other dimension

Types of Raster Data

  • Continuous: numeric data representing a measurement (e.g., elevation, precipitation)

  • Categorical: integer data representing factors (e.g., land use, land cover)

Continuous Rasters

mintemp <- rast("ftp://ftp.hafro.is/pub/data/rasters/Iceland_minbtemp.tif")
mintemp
class       : SpatRaster 
dimensions  : 340, 375, 1  (nrow, ncol, nlyr)
resolution  : 2000, 2000  (x, y)
extent      : -1050226, -300225.9, -699984.3, -19984.32  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=laea +lat_0=69 +lon_0=-4 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs 
source      : Iceland_minbtemp.tif 
name        : Iceland_minbtemp 
min value   :       -0.9982879 
max value   :        8.6031137

Categorical Rasters

  • Requires a classification matrix and coercion to factor

  • levels allows you to define the categories

# Create classification matrix
cm <- matrix(c(
  -2, 2, 0,
  2, 4, 1,
  4, 10, 2), ncol = 3, byrow = TRUE)

# Create a raster with integers
temp_reclass <- classify(mintemp, cm)
tempcats <- c("cold", "mild", "warm")
levels(temp_reclass) <- tempcats
cats(temp_reclass)
[[1]]
  value category
1     0     cold
2     1     mild
3     2     warm

Adding Dimensions

  • When data are aligned in space and/or time, more efficient to process as ‘cubes’ or ‘stacks’

  • Bands of satellite imagery, multiple predictors, spatio-temporal data

A note about support

  • We talked briefly about the agr option in the st_sf() function

  • agr refers to the attribute-geometry-relationship which can be:

    • constant = applies to every point in the geometry (lines and polygons are just lots of points)
    • identity = a value unique to a geometry
    • aggregate = a single value that integrates data across the geometry

  • Support is the area to which an attribute applies.

  • Rasters can have point (attribute refers to the cell center) or cell (attribute refers to an area similar to the pixel) support

Rasters in R

Rasters in R

  • raster - the original workhorse package; built on sp, rgdal, and rgeos
    • RasterLayer, RasterStack, and RasterBrick classes
  • terra - relatively new; developed by the raster folks, but designed to be much faster
    • SpatRaster and SpatVector classes
  • stars - developed by sf package developers; tidyverse compatible; designed for spatio-temporal data
    • stars class
    • Crosswalk between raster and stars is available here
    • Only way to deal with rectilinear and curvilinear data

Rasters by Construction

mtx <- matrix(1:16, nrow=4)
mtx
     [,1] [,2] [,3] [,4]
[1,]    1    5    9   13
[2,]    2    6   10   14
[3,]    3    7   11   15
[4,]    4    8   12   16
rstr <- terra::rast(mtx)
rstr
class       : SpatRaster 
dimensions  : 4, 4, 1  (nrow, ncol, nlyr)
resolution  : 1, 1  (x, y)
extent      : 0, 4, 0, 4  (xmin, xmax, ymin, ymax)
coord. ref. :  
source      : memory 
name        : lyr.1 
min value   :     1 
max value   :    16

Note: you must have raster or terra loaded for plot() to work on Rast* objects; otherwise you get Error in as.double(y) : cannot coerce type 'S4' to vector of type 'double'

Rasters by Construction: Origin

  • Origin defines the location of the intersection of the x and y axes
r <- rast(xmin=-4, xmax = 9.5, ncols=10)
r[] <- runif(ncell(r))
origin(r)
[1] 0.05 0.00
r2 <- r
origin(r2) <- c(2,2)

Rasters by Construction: Resolution

  • Geometry is implicit; the spatial extent and number of rows and columns define the cell size
  • Resolution (res) defines the length and width of an individual pixel
r <- rast(xmin=-4, xmax = 9.5, 
          ncols=10)
res(r)
[1] 1.35 1.00
r2 <- rast(xmin=-4, xmax = 5, 
           ncols=10)
res(r2)
[1] 0.9 1.0
r <- rast(xmin=-4, xmax = 9.5, 
          res=c(0.5,0.5))
ncol(r)
[1] 27
r2 <- rast(xmin=-4, xmax = 9.5, 
           res=c(5,5))
ncol(r2)
[1] 3

Rasters from Files

  • Building rasters useful for templates
  • More common to read from files
r <- rast(system.file("ex/elev.tif", package="terra"))
r
class       : SpatRaster 
dimensions  : 90, 95, 1  (nrow, ncol, nlyr)
resolution  : 0.008333333, 0.008333333  (x, y)
extent      : 5.741667, 6.533333, 49.44167, 50.19167  (xmin, xmax, ymin, ymax)
coord. ref. : lon/lat WGS 84 (EPSG:4326) 
source      : elev.tif 
name        : elevation 
min value   :       141 
max value   :       547

Accessing Raster Attributes

  • syntax is different for terra compared to sf